#include<stdio.h>
#include<stdlib.h>
#include<math.h>
#include<pthread.h>
#include<gmp.h>
#include<time.h>

#define NUMTHR 8

pthread_t t[NUMTHR];
mpf_t arg[18];

void * gauss_legendre(void *);

int main(void){
	int i = 0, cont=0;
	time_t t1, t2;
	mpf_t precisao;
	
	FILE *fp;
	
	fp = fopen("pi_input.txt", "r");
	if(fp == NULL){
		printf("Erro na abertura do arquivo!\n");
		return -1;
	}
	/*inicializando as variáveis*/
	for(i=0; i< 18; i++){
		mpf_init(arg[i]); 
	}	
	(void) time(&t1);
	/*inicializando valores*/
	mpf_set_d(arg[0], 1.0); 
	mpf_set_d(arg[1], 0.707106781); 	
	mpf_set_d(arg[2], 0.25); 	
	mpf_set_d(arg[3], 1.0);
 	/*define a precisao*/
mpf_init_set_ui(precisao, 10);
mpf_pow_ui(precisao, precisao, 10000000); 
mpf_ui_div(precisao, 1, precisao);

	printf("Retorno do input file = %d\n", mpf_inp_str (arg[9], fp, 10)); 
	
	/*Impressão na tela*/
/*	printf("Printing initialized variables:\n");
	mpf_out_str (stdout, 10, 10, arg[0]); printf ("\n");
	mpf_out_str (stdout, 10, 10, arg[1]); printf ("\n");
	mpf_out_str (stdout, 10, 10, arg[2]); printf ("\n");
	mpf_out_str (stdout, 10, 10, arg[3]); printf ("\n\n");
	mpf_out_str (stdout, 10, 10, arg[9]); printf ("\n\n");
*/

	while(mpf_cmp (arg[10], precisao) <= 0 || cont < 159){	
		for(i = 0; i< NUMTHR; i++ ){
			    pthread_create(&(t[i]), NULL, gauss_legendre, arg);
			    pthread_join(t[i], NULL);
		}
		cont+=i;
		mpf_sub(arg[10], arg[8], arg[9]); 
		mpf_abs(arg[10], arg[10]);
	}
	mpf_add(arg[10], arg[4], arg[5]); //(an + bn)
	mpf_pow_ui (arg[13], arg[10], 2); //(an + bn)²
	mpf_mul_ui (arg[17], arg[6], 4); //(4*arg[6])
	mpf_div (arg[8], arg[13], arg[17]); // termo de cima dividido pelo termo de baixo e armazena em pi

	(void) time(&t2);
	printf("Resultado final: Pi = "); mpf_out_str (stdout, 10, 10, arg[8]); printf ("\n");
	printf("Numero de iteracoes = %d, Tempo de calculo = %ld segundos\n", cont-1, (int) t2-t1);


	return 0;
}

void * gauss_legendre(void * argument){
	mpf_t *arg = (mpf_t*) argument;

	/*an = (a0+b0)/2;*/
	mpf_add(arg[10], arg[0], arg[1]);
	mpf_div_ui(arg[14], arg[10], 2.0);
	mpf_set(arg[4], arg[14]);

	/*bn = sqrt(a0*b0);*/
	mpf_mul (arg[12], arg[0], arg[1]);
	mpf_sqrt(arg[15], arg[12]);
	mpf_set(arg[5], arg[15]);	

	/*tn = t0 - p0*(a0 - an)*(a0 - an);*/
	mpf_sub (arg[12], arg[0], arg[4]);	
	mpf_pow_ui (arg[16], arg[12], 2);
	mpf_mul (arg[12], arg[3], arg[16]);
	mpf_sub (arg[16], arg[2], arg[12]);	
	mpf_set(arg[6], arg[16]);	

	/*pn = 2 * p0;*/
	mpf_mul_ui (arg[7], arg[3], 2);
	
	/*pi = (an + bn)*(an + bn)/(4*tn);*/	
	mpf_add(arg[10], arg[4], arg[5]);
	mpf_pow_ui (arg[13], arg[10], 2);
	mpf_mul_ui (arg[17], arg[6], 4);
	mpf_div (arg[8], arg[13], arg[17]);

	mpf_set(arg[0], arg[4]);
	mpf_set(arg[1], arg[5]);
	mpf_set(arg[2], arg[6]);
	mpf_set(arg[3], arg[7]);

//	printf("Pi (thread)= "); mpf_out_str (stdout, 10, 10, arg[8]); printf ("\n");	
}
